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Abstract Active Galactic Nuclei (AGN), such as 
Seyfert galaxies, quasars, etc., show light variations 
in all wavelength bands, with various amplitude and 
in many time scales. The variations usually look er- 
ratic, not periodic nor purely random. Many of these 
objects also show lognormal flux distribution and RMS 
- flux relation and power law frequency distribution. 
So far, the lognormal flux distribution of black hole 
objects is only observational facts without satisfactory 
explanation about the physical mechanism producing 
such distribution in the accretion disk. One of the most 
promising models based on cellular automaton mech- 
anism has been successful in reproducing PSD (Power 
Spectral Density) of the observed objects but could not 
reproduce lognormal flux distribution. Such distribu- 
tion requires the existence of underlying multiplicative 
process while the existing SOC models are based on 
additive processes. A modified SOC model based on 
cellular automaton mechanism for producing lognor- 
mal flux distribution is presented in this paper. The 
idea is that the energy released in the avalanche and 
diffusion in the accretion disk is not entirely emitted 
instantaneously as in the original cellular automaton 
model. Some part of the energy is kept in the disk 
and thus increase its energy content so that the next 
avalanche will be in higher energy condition and will 
release more energy. The later an avalanche occurs, 
the more amount of energy is emitted to the observers. 
This can provide multiplicative effects to the fiux and 
produces lognormal flux distribution. 
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1 Introduction 



Black hole objects, in a wide range of mass from Active 
Galactic Nuclei (AGNs) to black hole X-ray binaries 
(BHBs), generally show light fluctuation which is not 
periodic, nor purely random. This kind of variability 
is found in all wavelengths and in a wide range of time 
scales. It is one of the keys for understanding the nature 
of these objects. 

Long term monitoring of AGN light curves still fail to 
define clearly the upper limit of variability time scale, 
constraining our ability to extract inf ormation abou t 
disk structure from variability studies ( HawkinduOOm . 
For instance, twenty eight years of rnonito ring program 
of about 1500 quasars bv iHawkina ( 19961 ) still did not 
reach the common upper limit of quasars variability 
time scale. However, analysis of quasars long term 
light curve in quasars rest frame using Fou r ier po wer 
spectrum yield a break frequency. iHawkinsI (|2007l ) in- 
terpreted it as the 11 year characteristic time scale of 
quasars. 

Several attempts have bee n done to rep roduce 
such variability, for example, iTerrelll ( 1972 ) intro- 
duced a shot noise model, which was developed fur- 
ther by several authors (e.g. lLehtolll989tlLochner et al. 
19911 ). This model could reproduce the observed 



light curve and power spectral density (PSD) of the 
AGNs by superposition of the shots that occur ran- 
domly. In addition, the shots have identical time 
profiles with exponential decays. However, there is 
no cle a.r physical explanation for the creation of the 
shots (JTakeuchi and Mineshigd Il997l ). Other mod- 
els of variability that have been proposed so far, in- 
clu de starburst, di sk instability and microlensing mod- 
els ( Hawkinal2004 and references therein). 



Mineshige et al.l ( 1994al lb[) have tried to simulate the 



fluctuation in the light curve of black hole objects under 
the notion of self organized criticality (SOC) by using 
a cellular automaton rule. SOC is a continuous criti- 
cal state reached spontaneously by a system receiving 
energy and/or matter from outside, independent of its 
initial condition. 

The concept of SOC was first introduced by Bak et 
al. (1998), which used sand pile system as an illustra- 
tion. When grains of sand are dropped continuously in 
a surface, it will accumulate and create a peak which 
is higher and higher. When the peak is high enough 
and the slope of some part of the pile reach a critical 
value, then an additional falling grain of sand will trig- 
ger another grain to fall one after the other and creating 
an avalanche. After an avalanche the slope of the pile 
will be shallower. Additional sand drop will make the 
surface of the pile steeper and steeper, finally reach a 
critical value and trigger another avalanche, and so on. 
The sand pile will always be close to the critical state 
which may trigger avalanches. If the energy release of 
this sand pile system is recorded in a time series, it will 
produce l//-like va riability in its PSD. 

In their model, iMineshige et al.l ( 1994a ) and Mi- 
neshige et al. (1994b) divided the two-dimensional disk 
in radial and azimuthal directions, turned into cells. 
The disk receives mass input from a randomly chosen 
outermost cell and finally loses mass through the inner- 
most ring due to accretion onto its central black hole. 
Similar to the sand pile model, in the accretion disk 
the mass is transfcred inward through avalanche and 
diffusion processes that may happen in any cell under 
some special conditions. These processes convert gravi- 
tational potential energy into heat, eventually radiated 
from the disk. By using this model, they could success- 
fully simulate the light fluctuation and its PSD. The 
PSD showed a 1// tendency with power law index sim- 
ilar to those observed from some black hole objects. 

However, later observations of some AGNs revealed 
that som e of thern show lognormal X-ray flux distribu- 
tion (e.g. lGaskellll2004l) and f ollow hnear RMS-flux re- 
lation fe.g. lUttlev et al.ll2005l ). In this case RMS is the 
root mean square of the flux and the flux means the av- 
erage flux for the chosen bin. These properties are also 
found in some galactic acc reting source, such as Cyg 
X-1 and SAX J1808.4 -3658 dUttlev and McHardvll2001 : 
Gleissner et~aLll2004l) . The SOC model introduced by 



Mineshige et al.l (Il994allbll cou ld not explain such phe 



nomena (see e.g. lGaskellll2004l ). because, multiplicative 
processes are necessary to yield lognormal flux distri- 
bution, whereas the underlying process in the original 
SOC model is additive. What kind of physical processes 
in the accretion disk can yield lognormal flux distribu- 
tion? Some mathematical models have been introduced 



(see e.g. lGaskellll2004l:IUttlev et al.ll2005l) . but without 
a clear explanation on the physical mechanisms. 

In this work, we introduce an idea with reasonable 
physical mechanism which has multiplicative effect, in 
an attempt to obtain a lognormal flux distribution. 
We assume that the accretion flow is dominated by 
advective cooling instead of radiative cooling. This 
mechanism was called advection-dominated accretion 



flow (ADAF) (Ilchimarulll977l: lAbramowicz et al.lll988 : 



Kato et al.ll2008h 



In the present work, we adopt the idea that the en- 
ergy released by the viscous process is not efficiently 
emitted locally which could act as a multiplicative ef- 
fect that yield log normal flux distribution. I t is im por- 
tant to note that iTakeuchi and Mineshige! (|l997() cre- 
ated a numerical simulation of hydrodynamical model 
of advection dominated accretion disk with critical be- 
havior. They aimed to reproduce l//-like X-ray fluctu- 
ation with a physical mechanism focusing on the Cyg 
X-1 type fluctuation. We, on the other hand, aimed to 
address the lognormal flux distribution problem that 
was claimed to be unable to b e reproduced by a cellu- 
lar automaton model (see e.g. lGaskellll2004[ ). 

The plan of this paper is as follows: We will de- 
scribe the modified cellular automaton model in the 
second section in which we will also briefly review the 
lognormal distribution and the original cellular automa- 
ton model. In section three, we present the simulation 
details. The results and analysis will be presented in 
section 4. Section five is devoted for discussion and 
conclusion. 



2 The Model 

2.1 Lognormal Distribution 

The original SOC model, due to its additive nature, will 
have a normal flux distribution. The normal distribu- 
tion can be describe by the following formula: 
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where jj, and ct^ are the mean of the distribution and 
the variance of the distribution, respectively. 

A lognormal distribution is similar to a normal dis- 
tribution, but the variable is replaced by a logarithmic 
form 
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where r is the threshold parameter representing the 
lower limit on x. 

The lognormal distributions, which can describe 
many phenomena such as economic growth, popula- 
tion statistics, grain size in sand, etc., are also found 
in the X-r ay flux distri bution of some black hole ob- 
jects (e.g. iGaskelll 120041 ) . The fundamental difference 



between both distributions is their underlying process, 
additive in normal and multiplicative in lognormal dis- 
tribution. 

2.2 Advection Process 

In this work, in order to create multiplicative mech- 
anism in accretion disks of AGNs, we adopt the idea 
of the existing theory of advection process in accretion 
disk, which is well known as ADAF. 

Theoretically, there are two regimes of ADAF (see 
Kato et al.ll2008l for details) . The first regime is better 
known as radiatively inefficient accretion flow (RIAF) 
in which the radiation energy from the viscous process 
is not efficiently emitted locally (as in the case of the 
standard model) but stored as gas entropy and trans- 
ported inward. Due to inefficient radiative cooling, the 
gas becomes so hot that high-energy emission, in the 
form of hard X-ray, can be produced. This regime is 
also associated with a very low mass accretion rate. 
Such model was first developed by Ichimaru ( 19771) and 
extens ivel y discussed, later, b y lNarayan and Yil (|l994 , 



1995allb ). Abramowicz et al.l (J1995 . 1996D and others 



Naravan et al.l (jl995a ) showed that such model can re- 
produce the spectra of Sgr A*. It can also reproduce 
the spectra of some black hole binaries in their low lu- 
minosity state (also known as low/hard state). 

The other regime of ADAF occurs when the mass 
accretion rate is moderately high, close to or even ex- 
ceeding the critical mass accretion rate. Merit = Le/c^ 
and thus also known as supercritical accretion flow. 
Th e model of this ADAF re gime was first developed 
by lAbramowicz et al.l (|l988t ) 



at smaller radii and some fall onto the black hole with 
the accreting material without being radiated away. Su- 
percritical accretion has been observed in Galactic mi- 
croquasars, e.g. GRS1915-fl05. Some ultraluminous 
X-ray sources (ULXs), and narrow- line Seyfert 1 galax- 
ies (NLS ls) also seem to exhibit supercritical a ccretion 
(see e.g. IWatarai et al.ll200ll : lKawaguchill2003[ ) . 



2.3 Modified SOC Model 



In the original model ( Mineshige et al.lll994aM bl) an in- 
put of mass m is supplied at a randomly chosen cell in 
the outermost ring. From each ring there is mass dif- 
fusion m' from a randomly chosen cell to its adjacent 
inner cell. Suppose that i is index in radial and j in az- 
imuthal direction, if Mij is the mass of cell {i,j) where 
the diffusion occur, the mass will change according to 
the following rule : 
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Next, the model search the unstable sites which satisfy 
the instability criterion Mjj > Mcrit{r), where Mcritir) 
is a given critical mass. If, Mi_j > Merit then avalanche 
will occur according to the following rule : 
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and known as the slim 
disk model. Unlike RIAF, in slim disk model the ad- 
vective cooling does not always exceed the radiative 
cooling, e specially when M is bar ely exceed the criti- 
cal value ([Abramowicz et al.lll988( ). The important key 
for understanding supercritical accretion is the photon 
trapping effect. At high mass accretion rate, the disk 
becomes so optically thick that large number of inter- 
actions between photons and matter delay the radia- 
tion of energy released by the viscous process, espe- 
ciall y from the part near the ce ntral region of the disk 
(see lOhshuga et al.ll2003L 120051 ). That is, the radiation 
energy is trapped in the material flow and advected in- 
ward. Some of this trapped radiation energy is released ^ o^ 2^ 



The diffusion and the avalanche processes will release 
gravitational energy and produce radiation whose lumi- 
nosity prop ortional to the loss of gravita t ional poten- 
tial energy. iKuniava. Mineshige h, Ivezid (J2005I ) made 
a modification of the cellular automaton rule to incor- 
porate mass input not only via a cell but also all of 
the cells in the outermost ring receive small amount of 
matter from outside. 

The aim of the modification is to search for a mech- 
anism which can agree better to AGN observational 
facts. Unlike the case for the binary stellar mass 
black hole which receive mass input from the compan- 
ion through Lagrangian point only, it is more logical 
to think that mass input to the AGN's accretion disk 
comes from many sources surrounding the disk. 

The original SOC model could not yield lognor- 
mal flux distribution because of its intrinsic additive 
processes through avalanche and diffusion. For each 
avalanche or diffusion process, the energy released is 
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where Ami is the total amount of mass transported 
from the ring i to z + 1 in a time step by avalanche and 
diffusion. The energy produced from an avalanche in 
a cell is simply added to those of the other cells in the 
same time step, and assumed that the energy release is 
totally radiated locally. 

It is reasonable to think that not all of the energy 
produced are radiated locally. In other words, some 
fraction of energy remains in the disk and increases the 
energy content of the disk. Such m echanism has been 
consi dered in ADAF model (e.g. lAbramowicz et al. 
19881) . However it was not applied in the framework 



of SOC model. 

In the other avalanche, the part of the stored energy 
will also be radiated along with the currently produced 
energy. This will create a multiplicative effect needed 
to yield a lognormal flux distribution. So, when a suc- 
ceeding avalanche occurred in the same time step, the 
energy emitted is larger because the energy content of 
the disk has increased due to the absorption of some 
energy from the previous avalanche. 

Assume that the energy produced from an avalanche 
or diffusion is 
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Let the fraction of the energy emitted from an avalanche 
or diffusion is a, then the emitted energy is 



Li = aEi. 

The rest, Sj will be stored in the disk. 



(7) 



(8) 



In the next avalanche, still in the same time step, the 
energy produced is £'2 , a part of it will be emitted, also 
a part of the stored energy from the previous avalanche 
will be emitted: 

L2=Li+aE2+aSi. (9) 

The energy stored in the disk will be: 

52 = (1 - a)Si + (1 - a)E2 = (1 - a){Si + E2). (10) 

This energy will influence the amount of energy emitted 
in the next avalanche. Then, in general, the energy 
emitted after the n-th avalanche, will be : 



Ln — Ln-l + a{Sn-l + En). 

And the stored energy will be: 

^n = (1 - a){Sn-l + En). 



(11) 



(12) 



A black hole has a very strong gravity, even light can 
not escape from its gravitational pull, so that eventually 
there will be a certain amount of energy produced in the 
avalanches which will be swallowed by the black hole. 
In this model it is assumed that the amount of energy 
that falls into the black hole is equal to the amount of 
energy stored in the disk after the last avalanche in the 
same time step. Such mechanism occur in the slim disk 
model, in which the energy released by avalanche, diffu- 
sion or viscous flow is stored as rad iation entropy an d 
transported inward with accretion (JKato et al.l 120081) . 
The balancing between the remaining energy in the disk 
and the energy absorption by black hole may sound ar- 
tificial, but it is necessary to keep energy balance in the 
disk so that the disk always in a stable critical state. 



3 Simulation 

We did a numerical simulation based on cellular au- 
tomaton mechanism which was p reviously done by 
Kuniava. Mineshige fc Ivezig (|2005r ). The main aim of 



the simulation is to test whether it is possible to gen- 
erate lognormal flux distribution. Therefore we chose 
a set of fixed values of the previous parameters, that 
is m, the amount of each avalanche m', the amount of 
diffusion and m" , the amount of random mass input in 
the outermost ring of the accretion disk. 

Each time the avalanche or diffusion occurs, the 
amount of energy produced are calculated using equa- 
tion (5). The accumulated energy to be emitted after a 
series of n avalanche in a time step is calculated using 
equation (11) and the accumulated energy stored in the 
disk is calculated using equation (12). 

In the simulation, we tried several values of a, 0.01 < 
a < 0.8. The system with a < 0.5 means most of 
the dissipated energy from each avalanche or diffusion 
are not directly emitted but is stored in the disk and 
to be emitted part by part in the avalanches or diffu- 
sions afterwards at smaller radii. This setting is similar 
to those of advection dominated accretion flow model 



(|Kato et al.ll2008D . 

We ran the simulation for 10^ iterations. The first 
few thousands calculations were not used in the analysis 
to avoid the possibility that the system has not reach 
the SOC state. The sample of the simulated light curve 
is shown in figure 1. It is roughly similar to the typical 
light variation of black hole objects, that is irregular, 
not periodic nor purely random. We then analyzed light 
curve resulted from the simulation to find the clues of 
its fiux distribution. 
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Fig. 1 Sample of the light curve with a — 0.1 
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Fig. 2 Power Spectral Density of the simulated light curve 
with a — 0.1 



4 Result 

First we checked whether the simulated light curve 
have power law frequency distribution by calculating 
its PSD. Figure 2 shows an example of PSD for a = 0.1 
which has power law distribution, and is in good agree- 
ment with what were commonly derived from the ob- 
served light curve of AGNs. The PSD does not change 
much with the variation of a, and all yield similar slope. 
Another tools which can be used for checking the 
conformi ty of this model with observation is structure 
function (jSimonetti et al.lll985l ). So far, structure func- 
tion ha s been one of th e best analysis tool for AGN vari- 
ability ( Hawkinal2004 ) by which we can derive some in- 
formation from models and compare to those from the 
observerd light curve even for under sampled data and 
uneven sampling. Figure 3 shows the structure function 
derived from the simulated light curve. It shows a simi- 
lar characteristics with those generally derived from the 
observations of AGNs. 
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Fig. 3 Structure functions of the simulated light curve for 
several values of a 
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Fig. 4 Flux distribution for advection domination with 
a = 0.1, solid line is the theoretical lognormal function from 
equation 2 while the solid triangle points are from simulation 



To analyze the flux distribution of the simulated light 
curve, we grouped the flux into predefined bins and 
count the number of data in each bin. The width of 
the bin must be carefully determined, since if the bin 
is too wide then it will make the number of points in 
the distribution graph becomes too small while too nar- 
row bins will make the wing of the distribution under 
sampled. 

The flux distribution is then fitted to the lognormal 
function described in equation (2) by adjusting the pa- 
rameters T, a and fi. The flux distribution could fit well 
with the model of lognormal distribution with certain 
set of parameters and is shown in figure 4. 

For comparison, we present the flux distribution of 
the original SOC model without ADAF in figure 5. 
Roughly, it may also look fit, but when we study it more 
carefully, the points position relative to the line, then it 
is evident that the points in the wing arc slightly shifted 
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Fig. 5 Flux distribution for original SOC model (no ad- 
vection) with a — 1, solid line is the theoretical lognormal 
function from equation 2 while the solid triangle points are 
from simulation 
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Fig. 6 RMS - Flux relation derived from the simulated 
light curve, the solid line is the best fitted straight line to 
the data points 

to the left, while near the peak, the points shifted to 
right relative to the line. This means the points will 
be better fitted to a skewed distribution. The skewness 
of the other case with 0.1 < a < 1 is less clear com- 
pared to the a = 1 case. This shows that SOC model 
with cellular automaton mechanism can yield lognor- 
mal flux distribution by incorporating advection in the 
energy dissipation processes. 

In the light curve of many AGNs, RMS - flux relation 
arc frequently found. Such tendency can be seen also 
in our simulation with a = 0.1 (see figure 6) but not 
clearly appeared in the other calculation with larger 
values of a. This result makes us more convinced that 
the mechanism proposed here is in the right direction. 



5 Discussion 

The characteristics of simulated light curve which ap- 
pear in PSD, structure function, flux distribution, and 
RMS - flux relation are similar with those derived from 
observation. The slope of the PSD from our simu- 
lation is about 1.9, this is in the range of the PSD 
slope of AGNs (see for exampl e Van der Klid Il995 : 
Starling et al.ll2004t luttlev et al.|[2002l ). 

The structure function can provide some sort of in- 
formation about the variability time scale. The struc- 
ture function yi elded goes along wi th the observation 
(see for example lde Vries et al.ll2003l ) except the satura- 
tion part which provided the hint on longest variability 
time scale. Wc can see saturation part in our simula- 
tion, but not in the data of 30 years monit oring data of 
quasars reported bv Ide Vries et al.l (|2003l ). This indi- 
cates that the longest variability time scale of quasars 
is more than 30 years. 

The flux distribution of simulated light curve is 
lognor mal which is agree with observation. iGaskelll 
( 20041 ) found a lognormal flux distribution in an ex- 
treme narrow-line Seyfert 1 Galaxy, and argued that 
the lognormal intensity distribution of AGN is against 
the argument that the AGN variations being driven by 
self organized criticality. That is correct for the old 
version of SOC model. We, however successfully made 
some modification of the energy emission mechanism by 
combining ADAF and cellular automaton mechanism 
which can yield lognormal distribution, while keeping 
the other characteristics goes along with the observa- 
tion. 

This result can revive the SOC model for accretion 
disk which begin to fade away and more elaborate work 
based on SOC model is necessary to develop a better 
model with a hope to finally simulate the real situation 
in accretion disks. The model developed here is still 
rough and applied to the accretion disk in a not very 
detail manner, it's just to show that SOC model can 
still yield lognormal distribution and RMS-flux relation. 
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